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Abstract 

In  order  to  streamline  the  product  development  process,  the  design  space  for  FMVSS  201 U  impact  performance  of  a 
steel  mechanical  energy  absorber  assembly  was  investigated  by  means  of  LS-DYNA  970  explicit  finite  element 
simulation  methods  in  conjunction  with  statistical  analytical  procedures.  A  sequence  of  response  surfaces,  based  on 
various  levels  of  design  parameters,  was  developed  and  used  to  determine  minimal  stopping  distance  for  which  it 
would  be  possible  to  achieve  acceptable  impact  attenuation  performance  under  various  impact  loading  conditions 
given  a  worst  case  assumption  of  rigid  vehicle  interior  body  panels.  A  model  was  also  developed,  based  on  the 
variation  of  deterministic  variables,  in  order  to  estimate  and  minimize,  by  means  of  a  robustness  analysis,  the  range 
of  deviation  of  product  response  that  would  be  expected  as  a  result  of  variability  in  manufacturing  and  installation 

processes. 


Introduction 

Various  motor  vehicles  have  been  constructed  with  relatively  rigid  body  panels  and  various  other 
types  of  interior  assembly  that  offer  opportunities  for  reduction  of  the  incidence  and  severity  of 
head  injury  by  means  of  appropriately  designed  energy  absorbing  devices.  The  general 
performance  of  one  such  device,  a  steel  assembly  utilizing  plastic  deformation  of  fins  and  a 
cover  sheet  to  reduce  impact  severity  when  mounted  on  rigid  body  panels,  was  previously 
investigated  via  simulation  studies  and  it  was  determined  that  this  type  of  absorber  might  be  a 
worthwhile  candidate  for  application  to  the  solution  of  some  of  these  types  of  problem  [1]. 

As  a  continuation  of  that  work,  the  successive  response  surface  method  (SRSM)  of  Stander  and 
Craig  [2]  was  here  used  to  determine  the  characteristics  of  a  design  that  would  allow 
minimization  of  the  crush  space  required  to  enable  acceptably  low  levels  of  impact  severity. 

Classical  factorial  design  techniques  were  used  to  investigate  impact  response  robustness  in  the 
vicinity  of  the  optimum  combination  of  design  variable  settings.  It  was  found  that  a  small 
change  in  the  prescribed  center  point  of  one  of  the  design  variables  would  enable  a  significant 
improvement  in  the  robustness  of  the  design. 
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Energy  Absorber  Design  and  Simulation 

The  energy  absorber  system,  which  was  simulated  using  LS-DYNA  explicit  finite  element 
analysis  code  [3],  is  illustrated  in  Figure  1  and  was  comprised  of  a  network  of  0.5  inch  wide  mild 
steel  fins  that  were  sandwiched  between  a  rigid  panel  and  a  mild  steel  surface  panel  shell.  Each 
fin  had  a  0.25  inch  tab  at  each  end  for  use  as  an  attachment  to  the  surface  and  rigid  body  panels. 
The  fin  tabs  were  connected  to  one  another  by  means  of  a  mild  steel  web  with  a  view  toward 
ensuring  proper  fin  alignment  during  manufacture,  assembly,  and  installation.  The  fin  /  web 
assembly  was  connected  to  the  cover  sheet  via  spot  welds  and  was  constrained  at  the  rigid  panel 
end  using  single  point  constraints.  For  the  purposes  of  this  study  absorber  impact  performance 
was  adjusted  by  varying  crush  space,  spacing  distance  between  fins,  and  the  thickness  of  fins, 
webs,  and  cover  sheet. 


Figure  1.  Metal  fin  impact  energy  absorber. 


Component  level  impact  test  simulations  were  performed  with  the  intent  of  modeling  what 
would  be  expected  to  occur  during  US  Federal  Motor  Vehicle  Safety  Standard  (FMVSS)  201U 
upper  interior  tests  [4].  All  impact  simulations  involved  the  use  of  a  validated  First  Technology 
Safety  Systems  free  motion  headform  finite  element  model.  A  free  motion  headform  is  a 
modified  anthropometric  test  device  (crash  dummy)  head  that  is  instrumented  with  a  tri-axial 
accelerometer  in  order  to  measure  acceleration  during  the  course  of  impact. 

The  acceleration  history  from  each  impact  was  used  to  calculate  a  quantity  defined  as  the  head 
injury  criterion  (HIC)  according  to  Equation  (1),  where  a(t)  is  defined  as  the  resultant 
acceleration  as  a  function  of  time;  tl  and  t2  are  any  two  points  in  time  during  the  impact 
separated  by  not  more  than  36  milliseconds.  HIC(d)  is  a  correlation  between  HIC  for  the  free 
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motion  headform  and  HIC  for  a  full  50th  percentile  dummy  and  is  calculated  according  to 
Equation  (2).  Lower  HIC  is  better,  FMVSS  201 U  requires  that  HIC(d)  be  less  than  1000. 

For  all  simulations,  the  free  motion  headform  impacted  the  energy  absorber  with  an  initial 

["  Va(T)dr 

t  r  r /^t  I 

IL 

HIC{d)  =  0.75446  {Hid) +166.4  (2) 

velocity  of  1 5  mph  and  the  velocity  vector  of  the  head  at  impact  was  at  a  20°  angle  with  the 
normal  to  the  cover  surface.  Figure  2  illustrates  the  angle  of  the  velocity  vectors  for  various 
nodes  on  the  face  of  the  free  motion  headform. 
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Figure  2.  Angle  of  approach. 


The  successive  response  surface  method  (SRSM)  [2]  was  used  to  determine  the  combination  of 
three  design  factors  -  crush  space,  fin  spacing,  and  fm/web/cover  shell  thickness  -  that  would 
minimize  the  amount  of  crush  space  required  in  order  to  maintain  HIC(d)  <  700,  in  other  words, 
to  solve  the  optimization  problem 

minimize  crush  space 

subject  to  the  constraint  that 

HIC(d)  <  700. 


The  algorithm  used  here  closely  followed  the  default  implementation  of  the  successive  response 
surface  method  (SRSM)  as  outlined  in  LS-OPT  [5]. 

A  sequence  of  linear  response  surfaces  was  used.  33  full  factorial  designs  were  defined  to  be  the 
basis  experiments  for  each  of  the  response  surfaces.  D-optimal  subsets  -  with  seven  combinations 
of  each  of  the  three  design  factors  in  each  subset  -  were  selected  and  used  to  define  the  set  of 
simulations  that  were  performed  for  the  response  surfaces. 
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After  each  set  of  finite  element  analyses  was  performed,  regression  coefficients  were  calculated 
based  on  the  simulation  results  so  that  HIC(d)  could  be  approximated  as  a  linear  function  of 
crush  space,  fin  spacing,  and  fm/web/surface  shell  thickness.  Subsequently,  the  optimal 
combination  of  factor  levels  for  each  surface  was  determined  by  means  of  a  simple  sorting 
algorithm  and  used  as  the  center  point  for  the  next  iteration. 

Once  the  new  center  point  values  for  the  objective  function  -  crush  space  -  and  the  design 
variables  -  crush  space,  spacing,  and  shell  thickness  were  within  1%  of  those  of  the  current 
iteration,  the  process  was  considered  to  have  converged. 

If  convergence  was  not  achieved  subsequent  to  a  particular  iteration,  the  ranges  for  each  of  the 
design  variables  in  the  next  iteration  were  calculated  according  to  the  algorithm  outlined  in  [2] 
and  [5],  and  the  process  was  repeated. 


Figure  3.  Convergence  of  crush  space  and  fin  spacing. 
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Figure  4.  Convergence  of  shell  thickness. 

Convergence  was  quantitatively  achieved  after  15  iterations.  The  convergence  of  surface  center 
point  results  is  summarized  in  Figures  3  and  4  and  in  Table  1.  It  is  evident  from  the  results  that 
the  iteration  process  had,  for  the  most  part,  converged  to  optimal  levels  after  about  ten  iterations 

The  minimum  value  of  the  objective  variable  -  crush  space  -  was  determined  to  be  about  0.80 
inch  for  HIC(d)  to  be  less  than  700.  A  finite  element  simulation  at  the  predicted  optimum 
conditions  produced  a  result  of  HIC(d)  =  699. 


Initial 

Optimum 

Crush  space  (inch) 

0.875 

0.8044 

Fin  spacing  (inch) 

0.875 

0.9446 

Fin  /  web  /  cover  shell  thickness  (incl 

0.0285 

0.02616 

HIC(d) 

737 

699 

Table  1.  Convergence  of  design  input  and  output  parameters. 


The  acceleration  -  time  history  exhibited  by  the  optimum  design  is  shown  in  Figure  5. 
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Figure  5.  Acceleration  -  time  history  for  optimum  design. 


Robustness  Analysis  and  Design  Improvement 

Once  an  optimum  set  of  design  parameters  was  found,  an  analysis  was  performed  in  order  to 
estimate  the  degree  to  which  HIC(d)  performance  of  the  device  would  be  affected  by  moderate 
variations  in  crush  space,  web  spacing,  and  shell  thickness.  Some  degree  of  variability  in  HIC(d) 
response  would  be  expected  to  result  from  typical  manufacturing,  assembly,  and  installation 
processes.  An  improved  understanding  of  what  might  result  as  a  result  of  this  variability  would 
be  expected  to  help  to  ensure  greater  robustness  in  later  product  design  iterations. 

Myers  and  Montgomery  [6]  suggested  a  technique  for  investigating  robustness  involving 
analysis  of  the  interactions  between  noise  and  input  variables.  They  noted  that  interaction  plots, 
based  on  classical  factorial  experimental  methodology,  can  be  used  to  estimate  the  effects  of 
input  variable  interactions  on  the  mean  value  and  on  the  variability  of  output  variables. 

A  classical  23  full  factorial  design  with  each  of  the  design  factors  set  to  levels  ±5%  of  optimum 
was  used  to  develop  interaction  plots.  A  full  factorial  design  was  used  because  a  half  fraction 
would  not  offer  sufficient  resolution  to  prevent  aliasing  between  first  order  interactions. 

A  cube  plot  for  the  factorial  design  is  shown  in  Figure  6. 
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Figure  6.  Results  for  factorial  design. 

Two  of  the  design  variable  interactions  were  found  to  be  significant  -  the  interaction  between 
crush  space  and  shell  thickness  and  the  interaction  between  fin  spacing  and  shell  thickness 
(Figure  7).  Of  these  two,  the  interaction  between  crush  space  and  shell  thickness  was  found  to 
have  a  much  more  pronounced  effect  on  HIC(d)  response  than  the  interaction  between  spacing 
and  thickness.  The  variability  of  HIC(d)  with  changes  in  crush  space  is  reduced  somewhat  at 
higher  levels  of  shell  thickness. 


Figure  7.  Interaction  plot  for  crush  space  and  shell  thickness. 
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The  crush  space  -  shell  thickness  interaction  is  caused  by  the  combination  of  two  effects  (Figure 
8).  A  comparison  of  the  deceleration  peaks  that  occur  at  about  8  milliseconds  indicates  that  the 
peaks  for  the  devices  with  lower  shell  thickness  were  more  severe.  These  peaks  correspond  to  the 
increase  of  stiffness  that  is  related  to  the  onset  of  almost  complete  deflection  of  the  energy 
absorber  fins.  The  lower  level  of  shell  thickness  reduces  the  amount  of  free  motion  headform 
kinetic  energy  transformed  to  absorber  device  internal  strain  energy  during  the  fin  buckling 
phase  of  the  impact  event.  At  the  lower  level  of  crush  space,  as  the  head  exhausts  the  space 
available  for  deflection,  impact  severity  becomes  more  pronounced  as  a  result  of  increased  mean 
deceleration  due  to  reduced  stopping  distance. 


Figure  8.  Interaction  between  crush  space  and  shell  thickness:  acceleration  vs.  time  histories. 

At  the  higher  level  of  shell  thickness,  the  absorber  transforms,  during  fin  buckling,  significantly 
more  of  the  kinetic  energy  of  the  free  motion  headform  into  internal  strain  energy  of  the  device. 
As  a  result  the  disparity  between  HIC(d)  observed  at  the  two  levels  of  crush  space  is  much  less 
for  the  higher  level  of  shell  thickness  than  the  disparity  observed  for  the  two  crush  space  levels  at 
the  lower  level  of  shell  thickness. 

An  inspection  of  the  0.8446  inch  crush  space  (optimum  value  +5%)  face  of  the  cube  plot  for  the 
factorial  design  (Figure  6)  revealed  a  possible  opportunity  for  design  improvement  relative  to  the 
requirement  that  HIC(d)  be  less  than  700.  It  was  assumed,  based  on  the  results  from  the  factorial 
design,  that  HIC(d)  would  decrease  monotonically  with  increasing  crush  space  and  would 
probably,  for  values  of  crush  space  that  were  relatively  close  to  0.8446  inch,  be  minimized  for 
values  of  spacing  and  shell  thickness  that  were  somewhere  between  the  extreme  values  that 
were  set  in  the  factorial  design.  A  linear  interpolation  technique  was  used  to  iteratively  estimate  a 
minimum  value  of  crush  space  for  which,  given  spacing  fixed  at  0.9918  inch  and  shell  thickness 
set  at  0.02485  inch,  HIC(d)  would  be  expected  to  be  held  to  values  less  than  700. 

Using  this  technique,  a  new  value  for  nominal  crush  space  less  5%  was  determined  to  be  0.8562 
inch.  A  new  23  full  factorial  experiment  was  performed;  the  results  appeared  to  show  that  given  a 
new  nominal  value  for  crush  space  of  0.8975  inch  -  an  increase  of  less  than  0.01  inch  over  the 
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Figure  9.  Factorial  results  for  improved  design. 


original  nominal  value  -  the  new  design  would  indeed  meet  the  criterion  that  HIC(d)  be  less  than 
700  (Figure  9)  for  nominal  values  ±  5%. 


Inspection  of  the  interaction  plots  for  the  new  design  (Figure  10)  reveals  that,  with  the  new 
design,  there  is  very  little  interaction  between  crush  space  and  shell  thickness  or  crush  space  and 
spacing.  There  is  still  some  interaction  between  spacing  and  shell  thickness.  Although  an 


an  incremental  improvement  in  robustness  might  be  possible  with  an  increase  in  thickness,  it  was 
not  pursued  any  further  since  the  new  design  already  meets  the  criterion  of  HIC(d)  less  than  700 
for  nominal  values  of  design  variables  ±  5%. 

Further  characterization  of  the  behavior  of  the  new  design  was  achieved  by  development  of  a 
response  metamodel.  Comparative  results  reported  by  Simpson,  Lin,  and  Chen  [7]  revealed  that 
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uniform  designs  would  be  a  good  choice  for  a  space  filling  design  and  that  Kriging  would  be  an 
appropriate  choice  for  the  model  used  to  generate  a  response  surface.  Uniform  designs  based  on 
centered  L2  discrepancy  that  were  developed  using  the  threshold  accepting  method  [8]  were 
acquired  at  http://www.math.hkbu.edu.hk/UniformDesign/ . 

Kriging  models  are  of  the  form 


y=ZAfA*)+z(x) 


*=1 


(3) 


The  first  term  in  the  expression  is  a  regression  model;  the  second  term  is  a  model  of  a  random 
process  with  mean  zero  and  covariance 

V(x„Xj  )  =  cr2R(xi,xj )  (4) 

where  cT  is  the  process  variance  and  R  is  the  correlation.  A  Gaussian  spatial  correlation  function 
of  the  form 


R  ( xt  -  Xj )  =  exp  [  -  G  ( x,  -  Xj  )2  ] 


(5) 


was  used  for  the  purposes  of  the  current  work.  Kriging  calculations  were  performed  by  means  of 
the  DACE  MATLAB  Kriging  Toolbox  which  can  be  found  at 

http://www2.imm.dtu.dk/~hbn/dace/ .  The  technical  report  and  manual  [9]  associated  with  this 
Toolbox  give  a  detailed  explanation  of  the  calculations  associated  with  development  and  use  of  a 
Kriging  model. 

A  zero  order  polynomial  -  in  other  words,  a  constant  -  is  often  the  regression  term  of  choice  for 
Kriging  metamodels.  Martin  and  Simpson  [10]  demonstrated  and  suggested  the  value  of  the  use, 
at  times,  of  higher  order  polynomials  for  improving  Kriging  metamodel  approximation  of 
deterministic  computer  models.  Metamodels  were,  for  the  purposes  of  the  work  here,  developed 
using  uniform  designs  and  with  regression  models  comprised  of  a  constant,  of  a  linear 
polynomial,  and  of  a  quadratic  polynomial. 

Following  Simpson,  et.  al.  [7],  two  measures  were  used  to  attain  relative  estimates  of  goodness 
of  fit  for  the  metamodels.  The  first  involved  calculation  of  the  root  mean  squared  error  (RMSE) 
between  the  Kriging  estimate  and  the  known,  independent  results  generated  earlier  for  the 
improved  design  via  the  23  factorial  design  with  a  center  point. 


RMSE  = 


(6) 


The  second  measure  involved  calculation  of  the  maximum  error  between  the  metamodel  and  the 
known  points. 


Maximum  Error  =  max  |y;  -  yt 


(7) 
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Comparisons  of  the  various  combinations  of  sample  size  and  regression  model  are  given  in  Table 
2  in  terms  of  RMSE  and  in  Table  3  in  terms  of  maximum  error.  Due  to  the  fact  that  even  a 
simple  quadratic  regression  model  contains  10  terms  it  was  clearly  not  possible  to  generate  a 
Kriging  metamodel  with  a  quadratic  regression  for  a  sample  size  of  9. 


Sample  size 

Constant 

First  Order 
Polynomial 

Quadratic 

Polynomial 

9 

25.84 

19.47 

17 

12.35 

16.00 

19.02 

30 

23.95 

15.04 

10.67 

Table  2.  RMSE  comparison  for  various  Kriging  metamodels. 


Sample  size 

Constant 

First  Order 
Polynomial 

Quadratic 

Polynomial 

9 

42.70 

29.52 

17 

23.84 

37.44 

34.92 

30 

53.39 

30.08 

18.66 

Table  3.  Maximum  error  comparison  for  various  Kriging  metamodels. 


Although  there  were  no  clear  trends,  the  metamodel  that  minimized  these  measures  and  seemed 
to  thus  best  fit  the  known  factorial  results  was  the  one  generated  with  a  quadratic  polynomial  and 
with  a  sample  size  of  30. 

Results  from  that  metamodel  are  represented  in  contour  plots  of  HIC(d)  vs.  crush  space  and 
thickness  for  three  different  levels  of  spacing  (Figures  11,  12,  and  13).  The  ranges  for  the 
independent  variables  are  0.95-1 .05.  A  value  of  1.00  corresponds  to  the  nominal  value  for  the 
improved  design;  0.95  and  1.05  correspond  to  values  that  are  nominal  less  5%  and  nominal  plus 
5%,  respectively. 
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Figure  11.  Contour  plot  for  improved  design,  spacing  =  0.95  (5%  lower  than  nominal  value). 


Figure  12.  Contour  plot  for  improved  design,  spacing  =  1.00  (nominal  value). 


HIC(d) 
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Figure  13.  Contour  plot  for  improved  design,  spacing  =  1.05  (5%  higher  than  nominal  value). 


Slight  bias  error  between  the  Kriging  metamodel  and  the  known  values  from  the  factorial  design 
caused  contours  that  showed  HIC(d)  in  excess  of  700  near  combinations  of  0.95  crush  space, 
0.95  spacing,  1.05  shell  thickness  (Figure  11)  and  0.95  crush  space,  1.05  spacing,  0.95  shell 
thickness  (Figure  13).  Other  than  these  moderate  discrepancies,  the  metamodel  appeared  to 
support  a  thesis  that  HIC(d)  for  the  improved  design  remained  below  700  for  the  nominal  design 
values  ±  5%. 


Conclusions 

It  is  possible  to  very  efficiently  optimize  an  energy  absorber  design  using  the  LS-DYNA  explicit 
finite  element  code  in  conjunction  with  the  successive  response  surface  method  algorithm.  Use 
of  classic  factorial  techniques  in  combination  with  Kriging  response  surfaces  can  guide 
improvement  of  product  robustness  and  offer  insight  into  the  nature  of  a  product  and  its 
performance  variability.  An  enlightened  combination  of  these  techniques  enables,  if  nothing  else, 
valuable  and  relatively  inexpensive  insight  into  the  feasibility  and  behavior  of  various  design 
concepts. 
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